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ABSTRACT 

Context. Radial velocity measurements can be used to constrain the dynamical state of a stellar cluster. However, for clusters with 
velocity dispersions smaller than a few kms~' the observed radial velocity distribution tends to be dominated by the orbital motions 
of binaries rather than the stellar motions through the potential well of the cluster. 

Aims. Our goal is to characterize the intrinsic velocity distribution of a cluster from a single epoch of radial velocity data even for a 
cluster with a velocity dispersion of a fraction of a km s _1 . 

Methods. We investigate a maximum likelihood procedure, which was pioneered separately by Odenkirchen et al. (2002) and Kleyna 
et al. (2002). Assuming a period, mass ratio, and eccentricity distribution for the binaries in the observed cluster this procedure fits a 
dynamical model describing the velocity distribution for the single stars and center of masses of the binaries, simultaneously with the 
radial velocities caused by binary orbital motions, using all the information available in the observed velocity distribution. We test the 
capability of this procedure to reproduce the velocity dispersion of an observed cluster, using radial velocity data of an open cluster 
and Monte Carlo simulations. 

Results. We find that the fits to the intrinsic velocity distribution depend only weakly on the binary properties assumed, so the 
uncertainty in the fitted parameters tends to be dominated by statistical uncertainties. Based on a large suite of Monte Carlo simulations 
we provide an estimate of how these statistical uncertainties vary with the velocity dispersion, binary fraction, and the number of 
observed stars, which can be used to estimate the sample size needed to reach a specific accuracy. Finally we test the method on the 
well-studied open cluster NGC 188, showing that it can successfully reproduce a velocity dispersion of only 0.5 kms -1 using a single 
epoch of the multi-epoch radial velocity data. 

Conclusions. If the binary period, mass ratio, and eccentricity distribution of the observed stars are roughly known, this procedure 
can be used to correct for the effect of binary orbital motions on an observed velocity distribution. This allows for the study of the 
dynamical state of a stellar cluster with a small velocity dispersion from a single epoch of radial velocity data. 

Key words. Methods: statistical - open clusters and associations: general - Stars: kinematics and dynamics 



1. Introduction 

Most stars form in groupings or clusters, which are bound whilst 
still embedded in the birth molecular cloud and occasionally (~ 
10%) even after the gas has been expelled ( |Lada & La da 2003). 
In such an environment gravitational interactions between stars 
could dynamically process the binaries in the cluster, as well 
as disturb any disks or planetary systems orbiting the stars (e.g. 
|Kroupa| [1995] |Adams et ah] [2556} |Parker et al.| [20TT) . To fully 
understand the formation and evolution of these stellar clusters, 
as well as their effect on individual stars and multiple systems, 
we have to study their dynamical properties. 

For a resolved cluster with extensive photometric obser- 
vations, we can measure the surface density of the cluster. 
Assuming a density distribution along the line of sight for the 
stars (e.g. assuming a spherical cluster) we can use this to calcu- 
late the potential energy of the cluster. Combined with the radial 
velocity dispersion measured for an unbiased sample of stars, 
the cluster's kinetic energy can be calculated. The ratio of the 
kinetic to potential energy can be used to determine whether a 
cluster is either expanding, is collapsing, or is in virial equilib- 
rium and hence whether the cluster will be stable on the dynam- 
ical timescale. 

If the radial velocities of a large number of stars is available a 
more sophisticated dynamical model based on the potential en- 
ergy and density profile of the cluster can be tested. In such a 



model the velocity distribution is not characterized by a single 
velocity dispersion, but by a velocity distribution which varies 
with mass and/or distance from the cluster centre and the ve- 
locity distribution could be non-Gaussian. In globular clusters, 
measurements of the radial dependence of the velocity disper- 
sion have been used to for example show consistency of the stel- 
lar dynamics with the King (1965) thermal equilibrium model 
( |Gunn & Griffin||1979|, to constrain the anisotropy of the clus- 
ter (e.g. van de Ven et al. 2006 Sollima et al. 2012[) and even 
to check for deviations from Newtonian gravity (e.g. Baumgardt 
[etaL] [20051 |Sollima et aT] [2012] ). These studies were aided by 
the relatively high velocity dispersions and the relatively low bi- 
nary fractions (e.g. Gunn & Griffin 1979 Davis et al. 2008| > of 



most globular clusters. 

These studies can not be simply extended to lower mass 
clusters because these clusters tend to have smaller velocity 
dispersions. For small velocity dispersion (i.e. less than a few 
kms -1 ) the radial velocities from binary orbital motions can 
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1 (e.g. 


Kouwenhoven & de Grijs 


2008 


2009 


Gieles et al. 


McConnachie & Cote 201 Ok The importance of these 



binary orbital motions depends on the original velocity disper- 
sion, the binary fraction and the period, mass ratio, and eccen- 
tricity distributions of the binaries in the cluster. These veloci- 
ties induced by binary orbital motions completely mask the very 
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Fig. 1. Probability distribution function (left) and cumulative 
distribution function (right) for radial velocity offsets from the 
mean velocity for velocities drawn from a Gaussian distribution 
with a dispersion of 1 kms -1 (dotted, eq [2J, the velocity offsets 
expected due to binary orbital motions (dashed, eqs[3]and|4]i, and 
the velocity distribution for a cluster with a velocity dispersion 
of 1 km s -1 and a 100% binary fraction (solid), which is the con- 
volution of the dotted and dashed distribution (eq[T]i. Note that 
the binary orbital motions have a significant chance of giving a 
velocity offset of only tens of m/s, leading to a very high central 
peak in the probability density function, which has been cut off 
in the left panel. 



small velocity dispersions (< 1 km s ) expected in open clusters 



(e.g. Geller et al. 


2008 2010) and possibly in ultra-faint dwarf 


spheroidals (e.g. 


McConnachie & Cote 


2010 1. Even the larger 



velocity dispersions of young massive c usters ( e.g.|Bosch et al. 
20091 [Gieles et al.| |20l0l |Cottaar et"aTj [20T21 |Henault-Brunet 



et al.||2012 i, local star-forming regions (e.g .|Tobin et al.||2 009), 



and some low-mass globular clusters (e.g. Odenkirchen et aL[ 
[2002) |Blecha et al.| [20041 |Sollima et IE] |2012| i can be inflated 
by binaries. 

In Figure [T] we have illustrated the effect of binary orbital 
motions on the probability density function and cumulative dis- 
tribution of a single epoch of radial velocity data. The dotted 
line traces a Gaussian with a width of 1 km s , representing the 
velocity distribution of single stars and the center of masses of 
binaries in the observed cluster. The dashed line shows the prob- 
ability distribution of radial velocity offsets induced by binary 
orbits under several assumptions about the binary properties dis- 
cussed in Sec. [3] The solid line finally shows the convolution 
of these distributions, representing the distribution from which 
the observed radial velocities will be drawn, assuming a binary 
fraction of 100%. The binary orbital motions both broaden the 
Gaussian distribution and add a high-velocity tail to the distri- 
bution. The goal in this paper is to present a procedure which 
can be used after measuring a single epoch of radial velocities 
(having a distribution similar to the solid line) to characterize 
the intrinsic velocity distribution (dotted), assuming the binary 
period, mass ratio, and eccentricity distributions, which corre- 
sponds to fixing the shape of the dashed line. This fitting can 
accurately reproduce the intrinsic velocity distribution, because 
the velocities added due to binary orbital motions will generally 
have a very different velocity distribution from the intrinsic ve- 
locity distribution. 

An often used method to deal with the binary contamination 
is to obtain multi-epoch observations, which can be used to dis- 



tinguish between single and binary stars (e.g. Geller et al. 2008 
[20T0l |Tobin et al.[ 0O9\ |Cottaar et al] |2012b . The single stars 
(in this context generally defined as stars for which no significant 
radial velocity variations were found) as well as any binaries for 
which the velocity of the center of mass has been measured, can 



then be used to determine the velocity distribution of the cluster. 
However even after a large number of observations, there might 
still be undetected binaries, which broaden the measured veloc- 
ity distribution. For example a circular solar equal-mass binary 
with a period of a 1000 years, still has a velocity amplitude of 
2 kms -1 . 

In absence of multi-epoch data or to correct for any remain- 
ing undetected binaries Monte-Carlo simulations can be used to 
estimate for which velocity dispersion the broadening due to an 
assumed binary population would lead to the observed width of 



the distribution (e.g. 


Gieles et al. 2010 McConnachie & Cote 


2010; Sollima et al.l 


2012b. These methods tend to focus on re- 



producing only the width of the distribution, rather than the full 
velocity distribution observed. 

We investigate here an alternative strategy, which was used 
by Odenkirchen et al. ( 2002 ) in their study of the globular cluster 
Palomar 5 and separately by |Kleyna et al.| ( |2002[ l in their study 
of the Draco dwarf spheroidal. If the binary period, mass ratio, 
and eccentricity distributions are roughly known for the type of 
cluster being studied, this can be used to simultaneously fit the 
intrinsic velocity distribution and the fraction of binaries. The 
high-velocity tail caused by the orbital motions of close binaries 
can be used to determine the binary fraction in the observed clus- 
ter, with the distribution of velocities in the tail acting as a test 
on the assumed binary properties. From the high-velocity tail the 
effect of these binaries can then be extrapolated to lower veloc- 
ities to quantify the degree with which the velocity distribution 
has been broadened by binaries. 

In Sec. [2] we describe the maximum-likelihood analysis we 
use to fit the intrinsic velocity distribution and the binary orbital 
velocity distribution. In Sec. [4] we discuss the accuracy of the 
method, computing the statistical uncertainty as a function of 
the velocity dispersion, binary fraction, and number of observed 
stars, as well as the systematic offsets that can be introduced 
by differences between the assumed and true binary properties. 
Here we will also provide a tool to estimate the sample size 
needed to measure the dynamical state or binary fraction of a 
cluster with a desired precision prior to the observations. We test 
our model on the old open cluster NGC 188 in Sec[5] Geller et al. 



( 2008] |2009[ >; |Geller & Mathieu| ([201 1 ) present multi-epoch ra- 
dial velocity data for a large number of stars in this cluster and 
found that the single stars in their data set had a velocity disper- 
sion of 0.49^^ kms -1 . We show that we can reproduce their 
measured velocity dispersion using only a single epoch of their 
data. We discuss the limitations of the method and a possible 
extension to correct for any undetected binaries in multi-epoch 
data in Sec.[6]and summarize our results in Sec. [7] 

2. Method 

We consider a single epoch of iV bs radial velocity observations 
Vobs, i with measurement uncertainties of cr t> Si i of stars with 
mass estimates m,-. Our goal is to fit a dynamical model to the 
binary-corrected velocity distribution of the stars in the cluster 
for a set of assumptions about the binary period, mass ratio, and 
eccentricity distributions. As free parameters in this fit we will 
use the binary fraction in the cluster as well as any parameters 
used in the dynamical model (e.g. mean velocity, velocity dis- 
persion). 

First we calculate the likelihood function, which gives the 
likelihood of reproducing the observations given the parameters. 
This likelihood function is unique for every star, as it depends on 
the measurement error, the mass estimate of the star, and possi- 
bly other stellar properties taken into account in the dynamical 
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model (Sec. [XT}. Then we use the maximum likelihood estima- 
tor to convert the likelihood to reproduce the observations as a 
function of the parameters into a probability distribution of the 
parameters given the observed radial velocities (Sec. |2.2] >. 

2.1. Incorporating binaries in the velocity distribution 

The observed velocities can be seen as randomly drawn from 
the dynamical model (v'd yn ) describing the intrinsic velocity dis- 
tribution. We will take the measurement uncertainty into account 
directly in this dynamical model. So for the single stars this ve- 
locity will match the observed velocity. For a random subset of 
stars (the binaries) we add an additional velocity due to binary 
orbital motions (vb; n ). The probability for each star to be a binary 
depends on the binary fraction of the cluster and the distribution 
of the added velocities depends on the assumed period, mass ra- 
tio, and eccentricity distribution as well as the mass of the star. 
To calculate the probability of observing a certain velocity (v Q b s ) 
for a star in a binary, we have to integrate over the probability 
of any combination of Vd yn and Vbi n , for which Vd yn = v b s - Vbin- 
This corresponds to taking the convolution of the distributions of 
Vdyn and Vbin- The full velocity distribution is given by the sum 
of the distributions for single and binary stars: 



X,(v bs, i) = (1 -/bin)-£dyn, i(v bs, 0+ 

X+oo 
-Cdyn, i(Vobs, i ~ ^bin) ' -£bin, i(Vbin) dv b i n , 



(1) 



which gives the likelihood of measuring a velocity v b s , ; as a 
function of the binary fraction (/bin) and any other free parame- 
ters in the dynamical model. The subscript i is added to the like- 
lihood function (X;) to emphasize that this likelihood function 
is different for every star due to its dependence on the measure- 
ment error and the mass of the star. We note that in the deriva- 
tion of equation[JJwe assumed that the single stars and the center 
of masses of the binaries have the same dynamical distribution. 
Mathematically this assumptions is unnecessary and for some 
dynamical models it might be useful to drop this assumption 
(e.g. dynamical models that include mass segregation, where the 
velocity distribution depends on the total mass of a binary). 

In equation [JJ we have not taken any triple or higher-order 
systems into account. For solar mass these higher-order systems 
are rare (Raghavan et al. |2010| l. However, for more massive pri- 
maries the fraction of triple and higher-order system is higher 
( Zinnecker & Yorke 2007; Raghavan et al. 20101. For the sys- 
tem to be stable the radial velocity of any star will still be domi- 
nated by a single (generally the closest) companion, which could 
be either a single star or a multiple system with a much smaller 
semi-major axis. This implies that the radial velocity offsets ex- 
pected for a star in triple and higher-order systems can to first 
order still be approximated by that of a binary. 

In this paper we will focus on measuring the velocity dis- 
persion of a cluster. Using the maximum likelihood estimator 
the standard deviation and mean of a set of observed velocities 
can be derived by maximizing the likelihood that the observa- 
tions were drawn from a Gaussian distribution (|Pryor & Meylan| 



1993). So in order to measure the mean velocity and velocity dis- 
persion corrected for orbital motions of the stars in a cluster, we 
will use a Gaussian model to fit the intrinsic velocity dispersion 



Xdyn, i(^dyn) 



1 



: exp 



aj + cr 2 c ) 



2\ 



(Vdyn-^c) 



(2) 



which gives the likelihood for a velocity Vdyn given a measure- 
ment error cr, and the mean velocity fi c and velocity dispersion 
cr c of the cluster. 

If the observed stars have not been drawn randomly from the 
cluster, the velocity dispersion of the observed stars is no longer 
representative for the velocity dispersion of the whole cluster. 
As most spectroscopic observations will focus on the brightest 
stars, there will generally be a deviation between the observed 
velocity dispersion of these stars and the velocity dispersion of 
the cluster, if the cluster is mass segregated. When the observed 
stars are not representative of the whole population, a more real- 
istic model could still be fitted. In this case the velocity disper- 
sion (cr c ) should no longer taken to be a constant for all stars, 
but the dependence of the velocity dispersion on the distance 
from the cluster center or the mass of the star should be mod- 
eled. Another possible extension would be to model rotation, by 
allowing fi c to vary across the cluster. In our efforts to calculate 
the accuracy with which the intrinsic velocity distribution can be 
retrieved (Sec. |4]l we will focus on the simple case with only a 
Gaussian intrinsic velocity distribution and hence on the accu- 
racy with which we can measure a single velocity dispersion. In 
Sec.[5]we will have to adjust the model for the intrinsic velocity 
distribution when fitting the radial velocity data from NGC 188 
to take into account the contamination of fore- and background 
stars in the observed sample. 

Given a set of assumptions about the binary period, mass ra- 
tio, and eccentricity distributions, we numerically calculate the 
distribution of radial velocity offsets due to binary orbital mo- 
tions. We start by randomly drawing a large number (~ 10 5 ) of 
binaries with period Pj, mass ratio qj and eccentricity ej from 
the assumed distributions. In addition to these we draw for every 
binary a random time tj since the last passing of periastron. To 
minimize the time of the computation, we have pre-computed a 
dense two-dimensional grid giving the relative velocity between 
the binary components for various eccentricities (0 < e; < 1) and 
phases (0 < tj/Pj < 1) for a period of 1 year and a semi-major 
axis of 1 AU. The velocity taken from this grid can be used to 
calculate the velocity of the observed star through 



9j 



-Vgrid(e ptjlPj), 



(3) 



where to, is the mass of the observed star, and v gr id is the velocity 
calculated for the closest match in eccentricity and phase from 
the table. 

In order to calculate the distribution of velocity offsets, we 
still have to take into account the projection of the binary mo- 
tions along the line of sight. Such a projection for a star with an 
absolute velocity offset v ; in three dimensions will give a uni- 
form probability between -vj and +v, to get a certain velocity 
offset along the line of sight. This flat distribution implies that in 
the calculation of the total probability to have a velocity offset 
due to binary motions of Vb; n , the contribution of all binaries with 
v j < I vt,i n | is zero, while all binaries with vj > |vbi n | contribute a 
probability of ^ . For the list of primary velocities obtained (eq. 
[3j this gives for the likelihood of binary-induced motions to pro- 
duce a radial velocity offset v'bi n 



Xbii 



(4) 



jlVj>\v\jm\ 



where we sum the J- for all v ; > |vbi n |- The velocity distribution 
obtained using equations [3] and [4] for the solar-type field star bi- 
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nary property distributions discussed above has been shown in 
Figure [T] as a dashed line for a primary mass of 1 M Q . 

In practice we calculate the likelihood velocity distribution 
for binary-induced motions only for a single primary mass for a 
dense grid of possible radial velocity offsets. If the mass of an 
observed star is known, the likelihood distribution can be eas- 
ily converted to this new primary mass by multiplying the grid 

1 * 

of radial velocity offsets with (p^) 3 and dividing the computed 
probabilities by the same amount (keeping the normalization in- 
tact). If the period, mass ratio, or eccentricity distribution are 
known to change as a function of mass or some other observable 
(e.g. distance from the cluster center), we have to recalculate the 
likelihood function for binary orbital motions (eqs. [3]and|4jl for 
every observed star. 



2.2. Fitting the parameters 

Using equation [T] we can combine the intrinsic velocity distri- 
bution predicted by the dynamical model (e.g. eq. |2]i with the 
distribution of velocities due to binary orbital motions (eqs. [3] 
and [4]i to calculate the likelihood of observing a given velocity 
as a function of the binary fraction and any parameters describ- 
ing the dynamical model of the intrinsic velocity distribution. 
In its most general form this can be written as .£,(v;|x), where 
V; is the observed velocity, x is the set of free parameters, and 
X.i is the likelihood function, which will be different for every 
star, because it depends on the measurement error of the obser- 
vation, the mass of the star, and any other observables included 
in the dynamical model. The total likelihood of reproducing all 
observed velocities is then given by £(v|x) = rji-CiCViM- For a 
large number of observed stars this number will generally get too 
small for computers to compute (even for a well fitting model), 
so in practice we calculate the logarithm of the likelihood given 
byln(-C(v|x)) = 2: i ln(X ! (Vi|x)). 

The likelihood distribution of reproducing the observation 
given the parameters (X(v|x)), which has been computed above, 
has to be converted to the probability distribution of the free 
parameters given the observations (P(x|v)). Without any prior 
information on the probability distribution of the parameters, 
Bayesian inference states that £(v|x) and P(x|v) have the same 
distribution, with the only difference being is that P(x|v) has 
been normalized 



P(xW) ■ 



£(v|x) 
X-C(v|x)dx' 



(5) 



Because the two distributions have the same shape, the set of 
parameters which maximize the likelihood of reproducing the 
observations will be the best-fit parameters. 

As the number of fitted parameters increases, calculating the 
normalization in equation |5]becomes increasingly computation- 
ally expensive. An efficient way to explore the probability dis- 
tribution of the parameters is through a Markov Chain Monte 
Carlo simulation. This procedure will take a constrained ran- 
dom walk through the parameter space, producing a chain of 
parameters, whose distribution follows the probability distribu- 
tion. So the distribution of the values of a single parameter in 
this chain, will provide the probability distribution for this para- 
mater, which will be automatically properly marginalized over 
all the other free parameters. 



6000 
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Fig. 2. The distribution of absolute velocities (i.e. not projected 
along the line of sight) due to binary orbital motions calcu- 
lated for the period, mass ratio, and eccentricity distributions 



of the solar-type field stars (see text; Raghavan et al. 2010 



Reggiani & Meyer, 201 1 ) calculated from equation|3Jin blue and 
the best-fit log-normal in green. The log-normal has a width of 
0.841og 10 (kms _1 ) and a mean of 0.08 log 10 (kms ). 



3. Binary properties 

Because there is a degeneracy between the effects of the period 
distribution, mass ratio distribution, eccentricity distribution and 
binary fraction on the observed single-epoch velocity distribu- 
tion, it is not possible to fully constrain the binary properties 
using only a single epoch of radial velocity data. So a set of as- 
sumptions about these properties will have to be made to model 
the velocity distribution due to the orbital motions as described 
in Sec. |2.1| Because the binary properties have been shown to 



significantly change with mass ( Burgasser et al. 2007 Raghavan 
etal.||2010) Chini et al. 2012[ ), the choice of a period, mass ratio, 



and eccentricity distribution will generally depend on the mass 
range of the observed stars. 

In many clusters the observed stars will tend to be solar- 
type due to either the evolution of the more massive stars in 
dwarf spheroidals or low-mass globular clusters or due to the 
lack of massive stars in low-mass local star-forming regions. So 
throughout this paper we will use the well observed binary prop- 
erties for the solar-type field stars. For the period distribution we 
use the log-normal distribution with a mean period of 5.03 and a 
dispersion of 2.28 in logio days, which was found by Raghavan 
et al. (2010 1. For the s econdary to primary mas s ratio (q) we 



use the power-law from Reggiani & Meyer 



(2011) of 



dN 



1 



■0.5 



for 0.1 < q < 1. By setting the maximum q to one we make the 
assumption that the observed (brightest) star in a binary will al- 
ways be the most massive one. Raghavan et al. (2010 1 found a 



lack of high eccentricities for binaries with low periods. We take 
this into account be drawing the eccentricities from a flat proba- 
bility distribution between and e max with 



= V95+tanh(0.61og 10 -^ 
2 days 



-1-7)], 



(6) 



as proposed in Parker & Goodwin ( 2009} . For these binary prop- 
erties the distribution of the sizes of the three-dimensional veloc- 
ity differences between the primary star and the center of mass 
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of the binaries (calculated with eq. [3j can be very well approx- 
imated by a log-normal with a mean of 0.081og 10 kms -1 and a 
dispersion of 0.841og 10 kms _1 (see Figure]^. 

An important caveat with using the binary properties of the 
field population is that the solar-type binary properties might 
not be universal. The field star binary population might simply 
be a superposition of the binary properties in a large range of 
different star-forming regions with different binary properties. 
Furthermore we expect dynamical evolution to have altered the 
initial binary properties in older clusters (e.g.|Marks etal. 201 1 



Despite these caveats we show in Sec. [5] that we can accu- 
rately describe the single-epoch velocity distribution of the dy- 
namically old open cluster NGC 188 using the binary proper- 
ties of the solar-type field stars. In Sec. 4.2 we will discuss the 
systematic offsets in the measured velocity dispersion and bi- 
nary fraction induced, when the binary properties in the observed 
cluster do not match these assumptions. 



4. Accuracy of the procedure 

Besides the statistical (random) uncertainty, which can be cal- 
culated through the Markov Chain Monte Carlo simulations, 
there might be an additional systematic error in the best-fit in- 
trinsic velocity distribution if the assumed binary properties do 
not accurately describe the binaries in the cluster. Here we study 
the behavior of both the statistical and systematic uncertainty 
through a large suite of Monte Carlo simulations, where we cre- 
ate fake sets of radial velocity data and fit these with the proce- 
dure described in Sec. [2] 

In the Monte Carlo simulations described in Sec. l4.ll we en- 
sure that there is no systematic error by using the same binary 
properties to create the radial velocity data sets as we use to fit 
them, so we can focus purely on the statistical uncertainty and its 
dependence on the binary fraction, velocity dispersion, and num- 
ber of observed stars. In Sec. 4.2 on the other hand we greatly 



increase the number of observed stars to minimize the statistical 
uncertainty, so that we can focus on the effect on the measured 
velocity dispersion if we use different binary properties to create 
the data sets as we use to fit them. 

Throughout this section we will ignore the effect of the mea- 
surement errors of the observed radial velocities. For a normally 
distributed measurement error and a Gaussian intrinsic velocity 
distribution, the measurement error and the velocity dispersion 
are effectively indistinguishable (eq.pl. As long as the measure- 
ment error is sufficiently small that ~ cr? + °f> the effect of 
the measurement error on the observed velocity distribution is 
indeed negligible. For a larger measurement error we can still 
use the results from this section, however the velocity disper- 
sion should be replaced with the quadratic sum of the velocity 
dispersion and the measurement error. In this case the measure- 
ment error will have to be well quantified in order to retrieve the 
velocity dispersion. 

4.1. Statistical uncertainty 

After the velocities have been measured the statistical uncer- 
tainty can always be calculated through Markov Chain Monte 
Carlo simulations. When preparing for observations it is often 
useful to have an estimate of the accuracies that can be reached 
based only on the number of stars that will be observed and the 
expected properties of the cluster prior to any radial velocity 
measurements. 



With this goal in mind we ran Monte Carlo simulations on 
a large number of fake sets of radial velocities with solar type 
binary properties (Sec. [3]l with varying binary fractions, veloc- 
ity dispersions, and sample sizes. These data sets were created 
for stars with random primary masses between 0. 1 and 1 M Q 
following the Chabrier (200 5j) IM F, which is an update of the 
IMF presented in Chabrier (2003 ). For every star we calculated 
a random radial velocity from a normal distribution with the cho- 
sen velocity dispersion (eq. |2j. Every star has a probability set 
by the binary fraction to be part of a binary system, in which 
case an additional velocity due to the binary orbital motion is 
added. The binary properties are randomly chosen from the pe- 
riod, mass ratio, eccentricity and phase distributions and the bi- 
nary is randomly orientated with respect to the line of sight. The 
radial velocity dataset created in this way is then fitted using 
the procedure described in Sec. [2] Because we used the same 
assumptions about the intrinsic velocity distribution and the bi- 
nary properties to create and to fit the radial velocities, there are 
no systematic offsets in the measured velocity dispersion, mean 
velocity, and binary fraction with respect to the input velocity 
dispersion, mean velocity, and binary fraction. For every simu- 
lation the uncertainties on the best-fit parameters were calculated 
through Markov Chain Monte Carlo simulations. 

To first order (i.e. ignoring binaries) the uncertainty in 
the mean velocity is given by the standard deviation of the 



mean<5// c = 



fifes' 



where <x c is the velocity dispersion and N bs 



is the sample size. Similarly we expect the velocity dispersion 
to go as 6cr c == and finally the uncertainty in the binary 

fraction to go with 6 Am 7tU= ■ 

V^obs 

Using these first-order effects we can write the uncertainty 
in the binary fraction, velocity dispersion, and mean velocity as 



8 Am = 



do- c = 



(£/bin)c 



(£o"c)coit cr c 



(<5jUc)con- O"o 



(7a) 
(7b) 
(7c) 



obs 



where the corrections for the second-order effects due to binaries 
are given by (5/bin)con-, (£cr c ) corr , and (<5)u c ) corr . Figure [3] shows 
the dependence of these quantities on the velocity dispersion of 
the cluster (for A^bs = 300) and the number of observed stars (for 
cr c = 1 kms -1 ) for a binary fraction of 0% (blue, long dashes), 
50% (green, solid line), and 100% (red, short dashes). 

The bottom panels in Figure [3] show on the number of ob- 
served stars (A^ bs) for a fixed velocity dispersion of 1 km s -1 . We 
find that the second-order corrections on the velocity dispersion 
and mean velocity becomes independent of A^bs for N bs > 50. 
This means that the accuracy of the velocity dispersion and 
mean velocity can be well approximated to go with 1 / -yWobs- 
This dependence becomes unreliable for N bs < 50, because the 
measurement errors themselves vary greatly for small N bs an d 
the uncertainty in the velocity dispersion becomes strongly non- 
Gaussian. 

The uncertainty in the binary fraction (bottom left panel) 
shows a similar behavior, as long as the binary fraction is not 
too low. For a binary fraction of 50% (green) and velocity dis- 
persion of 1 kms -1 the binary fraction becomes effectively un- 
constrained between and 100% for N b s < 100. For a binary 
fraction of 0% we find that (£/bi n )corr declines as the sample size 
increases, implying that the measurement error in the best-fit bi- 
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Fig. 3. From left to right the second-order corrections of the uncertainties of the binary fraction (5/bi n )con-> the velocity dispersion 
(£>0"c)con> and the mean velocity ((5yu c )corr- The top panels show the dependency of these corrections on the velocity dispersion of 
the cluster (<t c ) and the lower panels show their dependency the sample size (W t>s)- These corrections have been computed for 
the a binary fraction of 0% (blue), 50% (green), and 100% (red). They can be used to calculate the one-sigma uncertainties using 
equations ' 



7a 



7c This allows for the estimation of the accuracy with which the intrinsic velocity distribution can be measured before 
any radial velocity data has been taken, as long as an estimate for the velocity dispersion and binary fraction are available. The cyan 
lines in the left panels show the lines with 6fbi n = 0.34, which is the expected uncertainty if the binary fraction is unconstrained (i.e. 
a flat probability distribution between and 1). 



nary fraction decrease more steeply than 1 / -yWobs over the full 
range of tested N \, s . 

Bee ause the second-order corrections (t5/bm)corri (6o~ c )corri 
and (^Uc)con- are constant over a broad range of A^bs it is possible 
to calculate the accuracy expected for future measurements from 
the top panels in Figure[3]using eqs 7a]|7c not only for N bs - 300 
for which the figure was calculated, but for any A^bs > 50 given a 
prior estimate of the velocity dispersion and binary fraction. To 
fit more complicated dynamical models with more free param- 
eters more stars will of course need to be observed to reach the 
same accuracy. 

For smaller sample sizes these estimates should be taken 
with caution. For small N \, s the uncertainties themselves be- 
come statistical quantities, which might randomly be much 
larger (or smaller) than expected. Furthermore the probability 
distribution of the parameters become non-Gaussian for small 
sample sizes. 



Keeping in mind these caveats eqs [7a] - 7c and the top pan- 
els in Figure |3]can be used to estimate the sample size needed to 
reach a certain accuracy on the observed velocity dispersion. For 
example if we consider a local star-forming regions with an ex- 
pected velocity dispersion of about 2 km s -1 and solar-type field 
star binary properties with a binary fraction of roughly 50%, we 
find from the upper center panel in Figure [3] (6cr c ) con = 1.1. 
To reach an accuracy of 10% on the velocity dispersion (i.e. 

2 



The calculations above were done assuming that the mea- 
surement uncertainty is significantly smaller than the velocity 
dispersion, so that the measurement uncertainty has a negligible 
effect on the observed velocity dispersion (effectively replacing 
cr\ + o"? « cp- c in eq. pi. Actual observations will sometimes in- 
clude stars for whicnthe measurement uncertainty is compara- 
ble or even larger than the velocity dispersion. These velocities 
will have a broader likelihood function with the width set by the 
quadratic sum of the measurement uncertainty and the velocity 
dispersion (eq. [2]). As seen in Figure [3] and eqs. 7a - 7c the un- 
certainty in observed parameters increases for broader distribu- 
tions, so these measurements will contribute less to the accuracy 
with which the intrinsic velocity dispersion and binary fraction 
are measured. However, as the measurement uncertainties are 
included in the model, these observations can be safely kept in 
the dataset as long as the measurement uncertainty is well mea- 
sured. If the measurement uncertainty is only a rough estimate it 
might be better to remove the observations with the largest mea- 
surement uncertainties from the dataset, as an overestimate of 
the measurement error might cause the velocity dispersion to be 
systematically underestimated and vice versa. 



4.2. The systematic uncertainty due to inaccurate binary 
assumptions 



-57 s = 0.1), we find from eq. 



7b 



would need to be observed to reach the desired accuracy. 



To fit an observed velocity distribution, we have to model the 
that A^ bs = ( qtj- J = 121 stars binary orbital motions. This requires assumptions for the binary 

period, mass ratio, and eccentricity distributions. If our assumed 
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Fig. 4. The ratio between the measured and input velocity dis- 
persion resulting from fits of artificial datasets with varying in- 
put velocity dispersions. These datasets were fitted assuming the 
binary properties of solar-type field star, but were created us- 
ing different binary period and mass ratio distributions. For the 
blue '+' -signs the binary population used to create the velocity 
distribution follows the same field solar-type binary properties 
used to fit the distribution, except that the mass ratio distribu- 
tion is flat rather than falling (see Sec. 4.2.1 1. For the green dots 



the massive stars follow the solar-type binary properties, while 
the binary properties of the low mass stars are set by the most 
likely binary fraction, period distribution, and mass ratio distri- 



bution observed for very low mass binaries from Burgasser et al. 
d2007l > (see Sec. |432] >. 



distributions do not match the actual binary properties in the ob- 
served cluster, there will be a systematic error in the measured 
velocity dispersion. 

Here we investigate the size of this systematic error by run- 
ning a set of Monte Carlo simulations. In these simulations we 
create artificial radial velocity datasets for a given input veloc- 
ity dispersion and a set of binary properties, which are different 
than those described for the solar-type field stars discussed in 
Sec. - 



2.1 



We fit these artificial radial velocity data sets assuming 
the binary period, mass ratio, and eccentricity distribution from 
the solar-type field stars discussed in Sec. 2.1 We then compare 



the measured velocity dispersion with the input velocity disper- 
sion used to build the artificial dataset and check for a systematic 
offset between the two. 

The measured velocity dispersion will differ from the input 
velocity dispersion with a random offset due to the statistical 
noise discussed in Sec. 4. 1 as well as a systematic offset due to 
the discrepancy between the binary properties used to create the 
artificial dataset and the binary properties used to fit the dataset. 
Because we are interested here in studying the systematic offset, 
we minimize the statistical noise by building a very large dataset 
of 10 5 radial velocities. Even with 10 5 observed stars there will 
still be some random scatter in the measured velocity dispersion 
and binary fraction. 

The procedure to build the artificial dataset is the same as in 
Sec. - 



4.1 



We briefly summarize it again. We draw the masses of 
the 'observed' stars following the Chabrier (2005 ) IMF between 
0.1 and 1 M . For every star a radial velocity is then drawn from 



a Gaussian distribution with a given input velocity dispersion. 
For the randomly assigned binaries a velocity offset due to the 
binary orbital motions is added. 

4.2.1. Flat mass ratio distribution 

To characterize the distribution of binary or bital motions we as 
sumed a mass ratio distribution of ^ ~ g~ 05 ( |Reggiani & Meyer 



|201 1 1. However, Reggiani & Meyer (20TT) found an uncertainty 
on this power-law slope of W ~ gr _u ;,±u - 2y . This implies that the 
observed mass ratio distribution is still consistent with being flat 
at the two sigma level. If the mass ratio distribution were in fact 
flat, this would cause a systematic offset in the measured veloc- 
ity dispersion given our assumption of a power-law mass ratio 
distribution. 

To explore this effect we ran a set of Monte Carlo simu- 
lations, where we create an artificial 'observed' radial velocity 
dataset assuming a flat mass ratio distributions between q — 0.1 
and q - 1 for several input velocity dispersions ranging from 0.1 
to lOkms -1 . To create these datasets we use the period and ec- 
centricity distributions of the solar-type field stars and fix the 
binary fraction to the solar-type field star binary fraction (46%: 
Raghavan et al. 2010| l. These artificial datasets are then fitted as- 
suming the (correct) solar-type field star period and eccentricity 
distributions, but with the (here incorrect) power-law slope for 
the mass ratio distribution. The velocity dispersion measured in 
this way have been plotted as blue '+'-signs in Figure|4]for val- 
ues of the input velocity dispersion between 0.1 and lOkms -1 . 

The flat mass ratio distribution has a higher fraction of near 
equal-mass binaries, so every individual binary has a larger prob- 
ability of creating a significant velocity offset than for the power 
law mass ratio distribution. So to reproduce the number of high 
velocity outliers in the artificial data set created with the flat mass 
ratio distribution, the fitted binary fraction will be higher than 
46% for the power-law mass ratio distribution. We find that the 
best-fit binary fraction is indeed higher than 46% and increases 
from 48% for an input velocity dispersion of 0.1 kms -1 to 55% 
for an input velocity dispersion of 10km s -1 . 

The effect of the flat mass ratio distribution is much smaller 
on the measured velocity dispersion than on the binary fraction. 
Figure [4] shows the ratio of the measured and input velocity dis- 
persion as blue '+' -signs. If the assumed binary properties would 
perfectly match those of the observed cluster this would always 
be one (except for some random scatter due to the statistical un- 
certainty discussed in Sec. |4.1| >. We see that the difference be- 
tween the flat and power-law slope for the mass ratio distribu- 
tion causes an offset of <1% to almost 2%. This offset increases 
for lower input velocity dispersion, because the effect of binary 
orbital motions become more important for narrower velocity 
distributions. 



4.2.2. Mass-dependent binary properties 

We also consider the case, where the low-mass binaries (primary 
mass between 0. 1 and 0.2 M Q ) in our artificial cluster do not have 
the same properties as the solar-type field stars, but rather the 
period distribution, mass ratio distribution, and binary fraction 
found by Burgasser et al. (2007 ) for very low mass binaries. The 
binary properties for these stars are still highly uncertain, but 
Burgasser et al. (2007 ) found that the observed binary properties 
for binaries with primary masses lower than 0.1 M Q could be 
fitted by a log-normal period distribution with a mean at 7 AU 
and a width of 0.24 log(AU), a power-law mass ratio distribution 
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with ^ ~ q 4S (for q < 0.6) and a binary fraction of 22%. This 
binary fraction is much lower than that observed for solar-type 
stars in the field, although the effect per binary on the velocity 
distribution is much larger, because they tend to be on closer 
orbits and tend to be equal-mass. 

In our simulations the masses range from 0.1 to 1 M . 
Although recent results (Janson et al. |2012| l suggest a smooth 
transition between the binary properties of the very low mass 
and solar-type stars, we will use a sharp transition between these 
regimes for this numerical test. For the low mass stars (< 0.2 M Q , 
~34% of all stars) we use the semi-major axis log-normal distri- 
bution, the mass ratio power-law distribution (for q between 0.1 
and 1), and the binary fraction found by Burgasser et al. (2 001) . 
For the higher mass stars we use the binary properties observed 
for the solar-type field stars (Sec. [3j with a binary fraction of 
46%. We once again fit these radial velocities under the (inac- 
curate) assumption that all stars still have the same binary frac- 
tion and follow the solar-type field binary period, mass ratio, and 
eccentricity distributions. These simulations are repeated for a 
broad range of input velocity dispersions. 

The measured velocity dispersion is shown in Figure [4] as 
green dots. Although the assumed solar-type field binary prop- 
erties are inaccurate for a third of the 'observed' stars, the offset 
between the measured and input velocity dispersion is still at 
most 4%. 



These two examples (Sec. 4.2.1 and 4.2.2 1 are meant to il- 
lustrate the effect of incorrect binary assumptions on the mea- 
sured intrinsic velocity distribution. For reasonable differences 
between the true and assumed binary properties, the effect on 
the intrinsic velocity dispersion turns out to be fairly small. For 
the cases discussed above only a few percent. It would require 
observations of thousands of stars, before this systematic offset 
becomes significant compared to the statistical error. However, 
if the binary properties are very uncertain, it might still be useful 
to check how much the intrinsic velocity distribution might vary 
for several reasonable assumptions about the binary properties 
HOdenkirchen et al.|[2002) i. 



5. Case study: NGC 188 

As a test for the procedure described above, we will here con- 

1999) 



sider the case of the old (7 ± 0.5 Gyr; Sarajedini et al. 
open cluster NGC 188. The dynamical state of the cluster was 
extensively stud ied by |Geller et al.] ( [20081 [2009) and |Geller & 



Math ieu] ( |201 \\ . Based on 9166 radial velocity measurements of 
1 108 stars in the direction of this cluster, Gell er et al.| {2008 ) 
identified a large number of stars with significant radial veloc- 
ity variations. Based on the remaining (seemingly single) 640 
stars, they derive a velocity dispersion of 0.64 ±0.04 kms -1 . This 
value is still inflated by unidentified binaries. |Geller & Mathieu 
( |201 l| l found that these unidentified binaries would inflate an 
intrinsic velocity dispersion of 0.49^ ^ kms -1 to the observed 
0.64 ± 0.04 kms -1 . In this section we will show that we can re- 
produce their intrinsic radial velocity distribution using only a 
single epoch of their data. 

We estimate the masses of the observed stars to be given 
by the closest point on the 7 Gyr Padova isochrone (Marigo 



etaLl [20081 with E(B-V) 



the photometric fit of Sarajedini et al 



0.09 an d DM : 
1999). 



11.44 (following 
The masses of the 

observed stars estimated in this manner vary between 0.8 and 
1.1 Mq. The sample of observed stars contains a significant con- 
tamination from fore- and background stars. For these stars we 
do not expect the computed masses to be accurate. Fortunately 
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Fig. 5. Histograms of the three single-epoch radial velocity dis- 
tribution from NGC 188 we consider in the top left, top right and 
bottom left panels, with the probability distribution of the best 
fitting model overplotted as the solid line for a star of 1 M and 
a measurement error of 0.4 kms -1 . The intrinsic velocity distri- 
butions (eq[8]l are plotted as dashed lines. In the bottom right the 
histogram of the velocities of the single stars identified in |Geller| 
|et al.| ( [2008] > is plotted. The intrinsic velocity distributions (renor- 
malized to the smaller population) computed from the three data 
sets have been overplotted. All three of these intrinsic velocity 
distributions are consistent with the velocity distribution of the 
single stars from |Geller et al. ( 2002) according to a KS-test (see 
text). In all panels the fraction of every bin expected to consist of 
non-members (Platais et al.l[2003) have been made yellow. The 
color coding of the lines corresponds to the data set used to fit 
the model. The first observed epoch has been shown in blue, the 
second epoch in green, and the last epoch in red. 



in our model the masses only affect the binary orbital motions, 
which have a negligible effect on the very broad velocity distri- 
bution of the fore- and background stars. 

These fore- and background stars will affect the intrinsic ve- 
locity distribution. Following |Geller et al.| ( |2008) we model the 
radial velocity distribution of the field star population with a 
Gaussian. Choosing the same distribution to fit the field stars al- 
lows us to make direct comparison with the results from |GelIer] 
|et al.| ( |2008) . We can add this Gaussian to the model for the in- 
trinsic velocity distribution (eq. [2]) to get 

' ( (Vdyn-yUc) 2 ^ 



-£dyn, i(^dyn) ~fc 



^2n(af+(r 2 c ) 
1 



exp 



d-/c) 



2(of + o- 2 )J 

(Vdyn-W) 



^2tt(o- 2 + o- 2 ) 



exp 



2 



2(crf + crj)) 



(8) 



where / c is the probability of an observed stars to be a clus- 
ter member (i.e. the fraction of cluster members) and o-f and /if 
are respectively the velocity dispersion and mean velocity of the 
non-members. These three parameters will be fitted in addition 
to the binary fraction (/bin), the cluster velocity dispersion (o" c ), 
and the mean velocity (// c ) of the cluster. 

Because multiple epochs of data are available for nearly all 
of the 1 108 stars we are able to run the procedure on three single- 
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epoch radial velocity data sets with only minimal duplication 
between the data sets. The radial velocities and measurement 
errors were taken from |Geller et al.| ( |2009| l. 

The first three panels in Figure [5] show the histograms of the 
radial velocities in these three data sets with the velocity distribu- 
tion for the best-fit model overplotted as solid lines. In principle 
there is no single best-fit distribution, because the velocity dis- 
tribution depends on the mass of the star and the measurement 
error. The distributions plotted in Figure [5] are for a star of 1 M Q 
and a measurement error of 0.4 km s -1 , which is typical for the 
vast majority of observed stars. We note that the model has been 
fitted to the individual observed velocities and is independent of 
the bin size chosen in Figure [5] 

The velocity distribution of single stars identified in the 
multi -epoch data by |Geller et al.| ( [2008| l is shown in the bottom 
right panel in Figure |5j where we have overplotted the intrin- 
sic velocity distributions found for the first data set (blue), the 
second data set (green), and the third data set (red). We check 
whether the observed radial velocity distribution of the single 
stars is consistent with being drawn from the predicted intrin- 
sic velocity distributions from the single epoch datasets using 
the Kolmogrov-Smirnov test. The KS test computes the largest 
deviation between the measured and predicted cumulative dis- 
tributions and calculates the probability that this deviation or a 
larger one can be created by chance when drawing random sam- 
ples from the predicted distribution. If this probability is small, 
the hypothesis that the measured and predicted distributions are 
the same can be rejected. 

Before we can compare the observed single star velocity dis- 
tribution with the intrinsic velocity distribution predicted from 
the single-epoch data, we have to correct for a selection effect. 
Because priority was given during the observations to proper 
motion members, the fraction of cluster members in the subset 
of observed stars with at least three epoch of observations over a 
baseline of one year is higher than in all observed stars. Because 



Geller et al.l d2008 ) require stars to have been observed for at 



least three epochs over at least one year before classifying them 
as single stars, the fraction of members among the stars identi- 
fied as single is also higher. Based on the proper motions mem- 
bership probabilities ( |Platais et al.| |2003) > we find the fraction 
of members among stars classified as single to be 53%, signifi- 
cantly higher than the 44% among all stars. 

After correcting for this selection effect by increasing / c to 
53% the KS test finds that the observed velocity distribution of 
the seemingly single stars is consistent with being drawn from 
the predicted intrinsic velocity distribution. The probabilities of 
drawing a distribution with at least the same offset in the cumu- 
lative distribution by chance is 17% for the first data set, 29% 
for the second data set, and 2% for the last data set. We note 
that even in the unlikely case that the assumptions going into the 
model (i.e. Gaussian distribution, field solar-like binary proper- 
ties) are perfect descriptions of NGC 188, the distributions are 
still expected to be slightly different, because the observed dis- 
tribution of stars identified as single might still be inflated by 
unidentified binaries. 

In Figure [6j we show the probability distributions for the six 
parameters we fit in our model. For comparison the literature val- 



ues for these parameters fr om|Platais et al. (2003 ); Geller et al. 
( |2008| l; |Geller & Mathieul ( |2011| > have been overplotted. Taking 



the average of all membership probabilities for all observed stars 
( |Platais et al.| |2003| l, we find that 44% of observed stars should 
be members. This is fully consistent with the results we find for 
the three individual data sets (top left panel of Figure 151. 
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Fig. 6. Posterior probability distributions of the free parameters 
in the model, enclosing 68.2%, 95.4%, and 99.7% (i.e. 1, 2, and 3 
cr) of the total probability. The six free parameters are the binary 
fraction (/bin), the fraction of cluster members (/ c ), the cluster 
velocity dispersion (<x c ), the cluster mean velocity (fi c ), the field 
velocity dispersion (cr f ) and the field mean velocity (//f). The 
different colors refer to the three different single-epoch radial 
velocity distributions with the same color coding as in Figure [5] 
The literature values for these parameters, derived from a proper 
motion analysis (Platais et al. 2003 i and the multi-epoch radial 
velocity data ( |Geller et al.||2008| |Geller & Mathieul|201 1[ > have 
been plotted as black dots. If the error bars are given, they have 
been plotted. 



For all three datasets we find a binary fraction consistent with 
being 100% (Fig.|6jl. This high binary fraction should be inter- 
preted with great caution, as we are only sensitive to the close 
binaries in NGC 188, which affect the observed velocity distri- 
bution. If a larger fraction of all binaries in NGC 188 is close 
compared to the solar-type field binaries, the total binary frac- 
tion might be significantly lower. Such a deviation from the bi- 
nary properties of the solar-type field stars might be caused by 
dynamical evolution, which decreases the pe riod of the hard bi 



naries (P < 10 b days in the case of NGC 188, Geller & Mathieu 



201 1 1. Evidence for significant dynamical evolution in NGC 188 
was found by von Hippel & Sarajedini ( 1998) and Bonatto et al. 
fl2003] l. 

Irrespective of whether there are wide orbit counterparts to 
the identified spectroscopic binaries, we can compare the bi- 
nary fraction found to the spectroscopic binary fraction from the 
multi-epoch data. |Geller & Mathieu| ( |2011| > found that 29 ± 3% 
of observed stars in NGC 188 are binaries with periods smaller 
than 10 4 days and q > 0.2, which makes these type of binaries 
2.4 + 0.2 times as frequent in NGC 188 as in the field. Under our 
assumptions that the binaries in NGC 188 have the same proper- 
ties as those in the field, the total binary fraction would then also 
be 2.4 ±0.2 times as high. This would lead to the unphysical con- 
clusion that 110+ 10% of stars in NGC 188 are binaries, which 
is consistent with the very high binary fractions found from the 
single-epoch data. The estimate of a total binary fraction above 
the 100% implies that the ratio of spectroscopic binaries to wider 
binaries is higher in NGC 188 than in the field. 
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We finall y compare the results of our Gaussian fits. |Geller| 
et al. (2008 1 measured a mean velocity of -35.31 kms and 
a dispersion of 38.01 kms -1 for the non-member s. The cluster 
memb ers have a mean velocity of -42.43 kms -1 ( Geller et al. 



2008) and a velocity dispersion of 0.49 + „™ kms 



0.08 



(Geller & 



Matnieu||2011 1. All of the mean velocities and velocity disper- 



sions measured for our three single-epoch datasets, except pos- 
sibly the velocity dispersion of the field population for which 



Geller et al. (2008 1 do not give an uncertainty, are consistent 



with the literature values within one to two sigma. 



6. Discussion 

Above we have shown that the presented likelihood procedure 
can accurately describe the intrinsic velocity distribution of 
NGC 188 and find an intrinsic velocity dispersion consistent 
with that the 0.49+^^ kms" 



from multi-epoch data (Geller & 



Math ieu||201 1[ ). This dispersion is accurately found, even though 
NGC 188 has a large spectroscopic binary fraction and there is a 
large contamination of fore- and background stars. 

In our analysis of NGC 188 we have attempted to use as little 
prior information about the cluster as possible in order to illus- 
trate the accuracy of the procedure based only on the observed 
single-epoch radial velocity and the mass estimated of the stars. 
The accuracy of the procedure could of course be improved by 
taking constraints from other sources into account. For example 
one could use the spectroscopic binary properties which [Geller| 



& Mathieu1 ( |2011| l found for NGC 188 instead of the solar-type 
field binary properties. Furthermore one could set the probabil- 
ity for a star to be a cluster member f c to the membership prob- 



abilities derived from the proper motion analysis of |Platais et al. 
(|2003|l. Fina lly one could use a model for the Milky Way (e.g. 
Robin et al. 2003 ) to predict the velocity distribution of the fore- 
and background star instead of fitting a Guassian to this distri- 
bution. 

Although we have tested the procedure here on an open clus- 
ter, it should be applicable to a broad range of stellar systems 
where the velocity distribution is significantly inflated due bi- 
nary orbital motions, such as faint Milky Way companion dwarf 
galaxies, young massive clusters, local star-forming regions, and 
some low-mass globular clusters. As long as there is a dynam- 
ical model for the velocity distribution of single stars and the 
center of masses of binaries and a rough approximation of the 
binary properties, a single epoch of radial velocity data can be 
used to constrain the dynamical state of the cluster. If enough 
velocity data is available it can serve as a test on the assumed 
binary properties for the observed cluster. 

6.1. Extension to multiple epochs 

Up till now we have only discussed the application of this pro- 
cedure to single-epoch data, because that is were the gain of us- 
ing this procedure is by far the largest. However, the method can 
also be applied to multi-epoch data. In theory the ideal procedure 
would be to compute for every star the probability distribution of 
possible center of mass velocities based on the full ensemble of 
binary orbital solutions, which are consistent with the data. This 
can then be used to constrain the intrinsic velocity distribution. 
However, if there are only a small number of epochs per star, 
calculating the full ensemble of binary orbital solutions would 
be computationally very expensive. 

An alternative procedure would be to only use those stars, 
whose radial velocity variations are consistent with them being 



single (e.g. |Geller et al.||2008||2010] l. The observed velocity dis- 
tribution of these stars would still be affected by any undetected 
binaries. The velocity offsets due to these binaries can be mod- 
eled by drawing a large sample of binaries with random periods, 
mass ratios, eccentricities, phases, and orientations. Any bina- 
ries which would have been detected given the timing of the ob- 
servations and the measurement errors are removed from this 



sample (e.g. |Geller et al. 2010| > . The distribution of velocity 
offsets caused by the remaining binaries, which would not be 
detected in the multi-epoch observations, is then set to the distri- 
bution of velocities due to binary orbital motions (Xbin, i( v bin))- 
This distribution can be used as input to fit the velocity distri- 
bution corrected for any undetected binaries, using exactly the 
same procedure as described in Sec. [2] 

6.2. Limitations 

The power of this procedure lies in the difference between the 
Gaussian intrinsic velocity distribution assumed and the high ve- 
locity outliers caused by the binary orbital motions. This implies 
that if a dynamical model is fitted, where the intrinsic velocity 
distribution is more similar to the distribution caused by binary 
orbital motions (i.e. has stronger wings), the accuracy of the pro- 
cedure would decrease. An example would be a model, which 
includes the possibility of stars ejected at a high velocity. To 
distinguish these ejected stars from velocity outliers caused by 
binary orbital motions, multiple epochs of radial velocity data 
will still be required. 

We have a similar problem in our analysis of the radial 
velocity data of NGC 188 (Sec. [5]l where the nature of the 
high-velocity outliers is unclear. They could either be members, 
whose velocity is offset due to its orbital motion, or they could 
be non-members. Due to this confusion the accuracy with which 
the binary fraction could be measured, was much lower as ex- 
pected for such a large set of radial velocity data. From the top 
left panel of Figure [3] we find that (5/bi n )corr for a cluster with 
a binary fraction of roughly 100% and a velocity dispersion of 
0.5 kms -1 is about 0.6. For the N * 500 members in NGC 188 



equation 7a gives an accuracy of 3%. This is much lower than 
the actual precision with which we could measure the binary 
fraction of ~ 15%. 

In addition to distinguishing between velocity outliers due 
to binary orbital motions and those caused by ejections or non- 
members, multi -epoch data has the advantage that it can be used 
to constrain the binary properties of the cluster. In principle the 
distribution of the high-velocity tail in the observed velocity 
distribution puts some constraints on the binary properties. If 
enough radial velocity data is available, this might be sufficient 
to reject a given set of binary assumptions. However the degen- 
eracy between the effects of the period, mass ratio, and eccen- 
tricity distributions and the binary fraction on this high-velocity 
tail means that these distributions can not actually be fitted. 

7. Summary 

We investigate here a maximum likelihood method to simultane- 
ously fit an intrinsic velocity distribution with the velocity off- 
sets induced due to binary orbital motions to a single epoch of 
observed radial velocities. The procedure was pioneered seper- 
ately by |Odenkirchen et al.| l |2002) and |Kleyna et al.| ( |2002) , but 
has not been used since. We show it can be used on other types of 
clusters, where it enables a proper treatment of the effect of bi- 
nary orbital motions on the radial velocities observed in surveys 
of single-epoch high-resolution spectra. 
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This procedure requires assumptions about the period, mass 
ratio, and eccentricity distribution for the binaries in the ob- 
served stars, although we show that reasonable uncertainties in 
the assumed binary properties do not significantly affect the in- 
trinsic velocity distribution we find (Sec. 4.2 1. 



Assuming a binary period, mass ratio, and eccentricity dis- 
tribution found for solar- type stars in the field (Ragh avan et al.| 
2010} the precision with which the binary fraction, velocity dis- 
persion, and mean velocity can be measured as a function of the 
number of observed stars, the velocity dispersion, and the binary 
fraction can be determined using eqs. 7a -|7c]and Figure[5] These 
results can be used to get a prior estimate of the number of stars 
that should be observed to reach a given accuracy. Although we 
calculate the accuracy with which the velocity dispersion can 
be measured, we emphasize that the procedure can be trivially 
expanded to test complicated dynamical models, where for ex- 
ample the velocity dispersion and the mean velocity depend on 
the position of the observed star in the cluster and/or its mass. 

As a test case we apply this procedure to three different sets 
of single epoch radial velocity data from 1 108 stars in NGC 188 
( Geller et al.j [2009). Assuming the period, mass ratio, and ec- 
centricity distribution from the field solar-type stars we are able 
to reproduce the velocity dispersion of ~ 0.5 km s for a single 
epoch of data, consistent with the value Geller & Mathieu ( |201 l| l 
derived from many epochs of radial velocity data, even though 
about half of the observed stars are fore- and background stars. 

For an assumed binary period, mass ratio, and eccentric- 
ity distribution this maximum likelihood procedure can be used 
to characterize the intrinsic velocity distribution for any set of 
observed radial velocities, where the observed distribution has 
been significantly broadened by binary orbital motions. This 
will happen in almost any clusters where the velocity disper- 
sion is lower than several kms -1 such as open clusters, faint 
Milky Way companion dwarf galaxies, young massive clusters, 
local star-forming regions, and some low-mass globular clus- 
ters. Determining the intrinsic velocity distribution will allow us 
to constrain the dynamical state of these clusters from a single 
epoch of radial velocity data, which gives a useful test on our 
understanding of the current state and evolution of these type of 
clusters. 
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